!==============================================================================!
!==                    GENETIC ALGORITHM PARALLEL PROGRAM                    ==!
!==============================================================================!
!==                                                                          ==!
!==  INSTITUTION: UTFPR - Federal Technological University of Paraná         ==!
!==  MAKER: Msc. Jonas Joacir Radtke                                         ==!
!==  E-MAIL: jonas.radtke@gmail.com                                          ==!
!==  LATTES: http://lattes.cnpq.br/7202527344656915                          ==!
!==  URL: http://paginapessoal.utfpr.edu.br/jonas                            ==!
!==  SVN CHECKOUT: https://gappp.googlecode.com/svn/trunk/                   ==!
!==                                                                          ==!
!==============================================================================!
program gapp

   implicit none

   include 'mpif.h'

   character(100) :: folderin = "./input_files/"    !< folder the for input files
   character(100) :: folderout = "./output_files/"  !< folder the for output files
   character(100) :: sname  !< simulations name
   character(100) :: ffit   !< name of executable for fitness calculation

   integer :: iarq = 10 !< number of input file
   integer :: reload    !< upload backup data
   integer :: g         !< generation
   integer :: ng        !< number of generations
   integer :: ps        !< population size
   integer :: nu        !< number of unknowns
   integer :: tl        !< total lenght of unknowns
   integer :: elitism   !< 0=off; 1=on
   integer :: ind       !< number of individual

   real(8) :: pc        !< crossover probability
   real(8) :: pm        !< mutation probability
   real(8) :: pmin      !< minimum probability in crossover
   real(8) :: rnd       !< random number

   integer, dimension(:), allocatable :: lu        !< lenght of unknowns

   integer, dimension(:,:), allocatable :: pop     !< population of chromosomes
   integer, dimension(:,:), allocatable :: newpop  !< new population of chromosomes
   integer, dimension(:,:), allocatable :: z       !< decoded to integer chromosome

   real(8), dimension(:), allocatable :: fit       !< fitness of the individuals
   real(8), dimension(:), allocatable :: rmin      !< minimum real value of the unknowns
   real(8), dimension(:), allocatable :: rmax      !< maximum real value of the unknowns

   real(8), dimension(:,:), allocatable :: r       !< decoded to real unknowns

   integer :: mate1     !< number of the individual for crossover
   integer :: mate2     !< number of the individual for crossover

   real(8) :: tcpu      !< total CPU time
   real(8) :: tcpuo     !< Accumulated CPU time before interruption
   real(8) :: tcpu1     !< starts CPU time
   real(8) :: tcpu2     !< finishes CPU time

   integer :: clock
   integer :: seed(97)
   integer :: iproc     !< mpi: identification number of the process
   integer :: code      !< mpi: status code
   integer :: nproc     !< mpi: number of processes
   integer :: tag       !< mpi: message label
   integer :: comm      !< mpi: MPI communicator
   integer :: status(mpi_status_size)  !< mpi: vector with information and source tag

   integer :: i
   integer :: iaux

   call mpi_init(code)

   call cpu_time(tcpu1)

   comm = mpi_comm_world
   tag = 42

   if (code /= mpi_success) then
      write(*,*) " =====  Erro ao iniciar o programa MPI. Encerrando...  ====="
      call mpi_abort(comm, code, code)
   endif

   call mpi_comm_size(comm, nproc, code)

   call mpi_comm_rank(comm, iproc, code)

   call system_clock(count = clock)
   seed = clock
   call random_seed(put = seed)

   if (iproc == 0) write(*,"(//, a, /)") &
      "  ==========  GENETIC ALGORITHM PARALLEL PROGRAM  ==========  "

   call get_parameters(folderin, folderout, iproc, reload, iarq, sname, tcpuo, &
      ffit, ps, ng, elitism, nu, tl, pc, pm, pmin, lu, pop, newpop, z, rmin, &
      rmax, fit, r)

   if (iproc == 0) then

      if (reload == 0) then

         open(20, file = trim(folderout) // trim(sname) // "-statistics.txt")
         51 format("# generation      meanfitness          fittest")
         write(20,51)

         open(21, file = trim(folderout) // trim(sname) // "-history.txt")
         52 format("# generation  individual          fitness       parameters")
         write(21,52)

      else

         open(20, position = "append", file = trim(folderout) // trim(sname) &
            // "-statistics.txt")
         open(21, position = "append", file = trim(folderout) // trim(sname) &
            // "-history.txt")

      end if

      call get_initial_population(ps, tl, iarq, reload, g, pop)

      do i = 1, nproc-1
         call mpi_send(g, 1, mpi_integer, i, tag, comm, code)
         call mpi_send(pop, ps*tl, mpi_integer, i, tag, comm, code)
      end do

   else

      call mpi_recv(g, 1, mpi_integer, 0, tag, comm, status, code)
      call mpi_recv(pop, ps*tl, mpi_integer, 0, tag, comm, status, code)

   end if

   if (g > ng) then
      write(*,"(//, a, //)") " === Limit generations already been reached! ==="
      stop
   end if

   call get_integer_value(ps, nu, tl, lu, pop, z)

   call get_real_value(ps, nu, rmin, rmax, lu, z, r)

   if (iproc == 0) then
      write(*,"(/, a, i4, a, /)") "Processing the", g, "th generation..."
   end if

   do ind = 1, ps

      iaux = mod(ind, nproc - 1) + 1

      if (iproc == iaux) then

         call get_fitness(folderout, ffit, ind, nu, r(ind,:), fit(ind))

         call mpi_send(ind, 1, mpi_integer, 0, tag, comm, code)
         call mpi_send(fit(ind), 1, mpi_double_precision, 0, tag, comm, code)

      end if

   end do

   call mpi_barrier(comm, code)

   if (iproc == 0) then

      do i = 1, ps

         iaux = mod(i, nproc - 1) + 1

         call mpi_recv(ind, 1, mpi_integer, iaux, tag, comm, status, code)
         call mpi_recv(fit(ind), 1, mpi_double_precision, iaux, tag, comm, &
            status, code)

         write(*,"(a, i4, a, 1pe15.8)") "The performance of the", ind, &
            "th individual is ", fit(ind)

         write(21,"(2(i12), 20(2x, 1pe15.8))") g, ind, fit(ind), r(ind,:)

         call flush(21)

      end do

      write(20,"(i12, 3(2x, 1pe15.8))") g, sum(fit)/ps, maxval(fit)

      call flush(20)

      call cpu_time(tcpu2)
      tcpu = tcpuo + tcpu2 - tcpu1

      call save_backup(folderout, sname, tcpu, g, ffit, ps, ng, elitism, nu, &
         tl, pc, pm, pmin, lu, pop, rmin, rmax)

   end if

   call mpi_barrier(comm, code)

   do g = g+1, ng

      if (iproc == 0) then

         write(*,"(/, a, i4, a, /)") "Processing the", g, "th generation..."

         do i = 1, ps, 2

            call selection_roullete2(ps, pmin, fit, mate1)
            call selection_roullete2(ps, pmin, fit, mate2)

            call random_number(rnd)

            if (rnd <= pc) then
               call crossover_simple(mate1, mate2, i, ps, tl, pop, newpop)
            else
               call nocrossover(mate1, mate2, i, ps, tl, pop, newpop)
            end if

         end do

         call mutation(ps, tl, pm, newpop)

         if (elitism == 1) call get_fittest(ps, tl, fit, pop, newpop)

         call replace(ps, tl, newpop, pop)

         do i = 1, nproc-1
            call mpi_send(pop, ps*tl, mpi_integer, i, tag, comm, code)
         end do

      else

         call mpi_recv(pop, ps*tl, mpi_integer, 0, tag, comm, status, code)

      end if

      call get_integer_value(ps, nu, tl, lu, pop, z)

      call get_real_value(ps, nu, rmin, rmax, lu, z, r)

      do ind = 1, ps

         iaux = mod(ind, nproc - 1) + 1

         if (iproc == iaux) then

            call get_fitness(folderout, ffit, ind, nu, r(ind,:), fit(ind))

            call mpi_send(ind, 1, mpi_integer, 0, tag, comm, code)
            call mpi_send(fit(ind), 1, mpi_double_precision, 0, tag, comm, code)

         end if

      end do

      call mpi_barrier(comm, code)

      if (iproc == 0) then

         do i = 1, ps

            iaux = mod(i, nproc - 1) + 1

            call mpi_recv(ind, 1, mpi_integer, iaux, tag, comm, status, code)
            call mpi_recv(fit(ind), 1, mpi_double_precision, iaux, tag, comm, &
               status, code)

            write(*,"(a, i4, a, 1pe15.8)") "The performance of the", ind, &
               "th individual is ", fit(ind)

            write(21,"(2(i12), 20(2x, 1pe15.8))") g, ind, fit(ind), r(ind,:)

            call flush(21)

         end do

         write(20,"(i12, 3(2x, 1pe15.8))") g, sum(fit)/ps, maxval(fit)

         call flush(20)

         call cpu_time(tcpu2)
         tcpu = tcpuo + tcpu2 - tcpu1

         call save_backup(folderout, sname, tcpu, g, ffit, ps, ng, elitism, &
            nu, tl, pc, pm, pmin, lu, pop, rmin, rmax)

      end if

      call mpi_barrier(comm, code)

   end do

   if (iproc == 0) then

      close(20)
      close(21)

      call get_fittest(ps, tl, fit, pop, newpop)

      call plot_statistics(sname, folderout)

      call plot_history(sname, folderout, nu, ps, ng, rmin, rmax)

      call cpu_time(tcpu2)
      tcpu = tcpuo + tcpu2 - tcpu1

      call write_results(tcpu, nu, ps, fit, r)

   end if

   call mpi_finalize(code)

contains

   !============================================================================

   subroutine get_parameters(folderin, folderout, iproc, reload, iarq, sname, &
         tcpuo, ffit, ps, ng, elitism, nu, tl, pc, pm, pmin, lu, pop, newpop, &
         z, rmin, rmax, fit, r)
      implicit none
      character(100), intent(in) :: folderin  !< folder for the input files
      character(100), intent(in) :: folderout !< folder for the output files
      character(100), intent(out) :: sname    !< simulations name
      character(100), intent(out) :: ffit     !< name of executable for fitness calculation

      integer, intent(in) :: iproc     !< identification number of the process
      integer, intent(inout) :: iarq   !< number of input file
      integer, intent(out) :: reload   !< upload backup data
      integer, intent(out) :: ps       !< population size
      integer, intent(out) :: ng       !< number of generations
      integer, intent(out) :: elitism  !< 0=off; 1=on
      integer, intent(out) :: nu       !< number of unknowns
      integer, intent(out) :: tl       !< total lenght of unknowns

      real(8), intent(out) :: tcpuo    !< accumulated CPU time before interruption
      real(8), intent(out) :: pc       !< crossover probability
      real(8), intent(out) :: pm       !< mutation probability
      real(8), intent(out) :: pmin     !< minimum probability in crossover

      integer, dimension(:), allocatable, intent(out) :: lu       !< lenght of unknowns
      integer, dimension(:,:), allocatable, intent(out) :: pop    !< population of chromosomes
      integer, dimension(:,:), allocatable, intent(out) :: newpop !< new population of chromosomes
      integer, dimension(:,:), allocatable, intent(out) :: z      !< decoded to integer chromosome

      real(8), dimension(:), allocatable, intent(out) :: rmin     !< minimum real value of the unknowns
      real(8), dimension(:), allocatable, intent(out) :: rmax     !< maximum real value of the unknowns
      real(8), dimension(:), allocatable, intent(out) :: fit      !< fitness of chromosome

      real(8), dimension(:,:), allocatable, intent(out) :: r   !< decoded to real unknowns

      integer :: i

      open(iarq, file = trim(folderin) // "input_file.txt")

      read(iarq,*) sname
      read(iarq,*) reload

      if (reload /= 0) then
         iarq = 11
         open(iarq, file = trim(folderout) // trim(sname) // "-backup.txt")
         read(iarq,*) sname
         read(iarq,*) tcpuo
      else
         tcpuo = 0.d0
      end if

      read(iarq,*) ffit
      read(iarq,*) ps
      read(iarq,*) ng
      read(iarq,*) pc
      read(iarq,*) pm
      read(iarq,*) pmin
      read(iarq,*) elitism
      read(iarq,*) nu

      sname = adjustl(sname)
      ffit = adjustl(ffit)

      allocate(lu(nu))
      allocate(rmin(nu))
      allocate(rmax(nu))

      tl = 0
      do i = 1, nu
         read(iarq,*) lu(i)
         read(iarq,*) rmin(i)
         read(iarq,*) rmax(i)
         tl = tl + lu(i)
      end do

      close(10)

      ! pm = pm/tl

      allocate(pop(ps,tl))
      allocate(newpop(ps,tl))
      allocate(z(ps,nu))
      allocate(fit(ps))
      allocate(r(ps,nu))

      pop = 0
      newpop = 0
      z = 0
      fit = 0.d0
      r = 0.d0

      if (iproc == 0) then
         open(22, file = trim(folderout) // trim(sname) // "-summary.txt")

         write(22,*)
         write(22,*) " =====  GENETIC ALGORITHM PARALLEL PROGRAM  ====="
         write(22,*)

         if (len(trim(sname)) <= 15) then
            write(22,"(a15, a)") trim(sname), " = sname:   Simulation name"
         else
            write(22,"(a, a)") trim(sname), " = sname: Simulation name"
         end if

         if (len(trim(ffit)) <= 15) then
            write(22,"(a15, a)") trim(ffit), &
               " = ffit:    Name of executable for fitness calculation"
         else
            write(22,"(a, a)") trim(ffit), &
               " = ffit: Name of executable for fitness calculation"
         end if

         write(22,"(i15, a)") ps, " = ps:      Population size (even number)"
         write(22,"(i15, a)") ng, " = ng:      Number of generations"
         write(22,"(1pe15.8, a)") pc, " = pc:      Crossover probability [0:1]"
         write(22,"(1pe15.8, a)") pm, &
            " = pm:      Mutation probability / (lu(1) + ... + lu(nu)) [0:1]"
         write(22,"(1pe15.8, a)") pmin, &
            " = pmin:    Minimum probability in crossover [0:1]"
         write(22,"(i15, a)") elitism, " = elitism: (0=off; 1=on)"
         write(22,"(i15, a)") nu, " = nu:      Number of unknowns"

         do i = 1, nu
            write(22,"(i15, a, i1, a, i1, a)") lu(i), &
               " = lu(", i, "):    Length of the ", i, " unknown"
            write(22,"(1pe15.8, a, i1, a, i1, a)") rmin(i), &
               " = rmin(", i, "): Lower value for the ", i, " unknown"
            write(22,"(1pe15.8, a, i1, a, i1, a)") rmax(i), &
               " = rmax(", i, "): Biggest value for the ", i, " unknown"
         end do
         write(22,*)
      end if

   end subroutine get_parameters

   !============================================================================

   subroutine get_initial_population(ps, tl, iarq, reload, g, pop) ! two last parameter are output
      implicit none
      integer, intent(in) :: ps           !< population size
      integer, intent(in) :: tl           !< total lenght of unknowns
      integer, intent(in) :: iarq         !< number of input file
      integer, intent(in) :: reload       !< upload backup data
      integer, intent(out) :: g           !< generation
      integer, intent(out) :: pop(ps,tl)  !< population of chromosomes

      integer :: i, j
      real(8) :: rnd

      if (reload == 0) then

         g = 1

         do i = 1, ps

            do j = 1, tl
               call random_number(rnd)
               pop(i,j) = int(2*rnd)
            end do

         end do

      else

         read(iarq,*) g
         g = g + 1

         do i = 1, ps
            read(iarq,*) pop(i,:)
         end do

         close(iarq)

      end if

   end subroutine get_initial_population

   !============================================================================

   !> \brief Get the fitness of each chromosome
   subroutine get_fitness(folderout, ffit, ind, nu, par, fit) ! last parameter is output
      implicit none
      character(100), intent(in) :: folderout !< folder the for output files
      character(100), intent(in) :: ffit      !< name of executable for fitness calculation
      integer, intent(in) :: ind    !< individual number
      integer, intent(in) :: nu     !< number of unknowns
      real(8), intent(in) :: par(nu)!< parameters for optimization
      real(8), intent(out) :: fit   !< fitness

      character*3 :: char3
      character(100) :: arqpar !< name of file for the parameters
      character(100) :: arqfit !< name of file for the fitness
      integer :: i

      call convert_int_to_char3(ind, char3)

      arqpar = trim(folderout) // "parameters" // char3 // ".txt"

      open(27, file = trim(arqpar))

      write(27,"(i15, a)") ind, " = ind: individual number"
      do i = 1, nu
         write(27,"(1pe15.8, a, i1, a)") par(i), " = par(", i, ")"
      end do

      close(27)

      call system("./" // trim(ffit) // " < " // trim(arqpar) // " > /dev/null")

      arqfit = trim(folderout) // "fitness" // char3 // ".txt"

      open(26, file = trim(arqfit))
      read(26,*) fit
      close(26)

   end subroutine get_fitness

   !============================================================================

   subroutine get_real_value(ps, nu, rmin, rmax, lu, z, r) ! last parameter is output
      implicit none
      integer, intent(in) :: ps        !< population size
      integer, intent(in) :: nu        !< number of unknowns
      integer, intent(in) :: lu(nu)    !< lenght of unknowns
      integer, intent(in) :: z(ps,nu)  !< decoded to integer chromosome
      real(8), intent(in) :: rmin(nu)  !< minimum real value of the unknowns
      real(8), intent(in) :: rmax(nu)  !< maximum real value of the unknowns
      real(8), intent(out) :: r(ps,nu) !< decoded to real unknowns

      integer :: i, j

      do i = 1, ps
         do j = 1, nu
            r(i,j) = (rmax(j) - rmin(j))/(2.d0**dble(lu(j)) &
               - 1.d0)*dble(z(i,j)) + rmin(j)
         end do
      end do

   end subroutine get_real_value

   !============================================================================

   subroutine get_integer_value(ps, nu, tl, lu, pop, z) ! last parameter is output
      implicit none
      integer, intent(in) :: ps           !< population size
      integer, intent(in) :: nu           !< number of unknowns
      integer, intent(in) :: tl           !< total lenght of unknowns
      integer, intent(in) :: lu(nu)       !< lenght of unknowns
      integer, intent(in) :: pop(ps,tl)   !< population of chromosomes
      integer, intent(out) :: z(ps,nu)    !< decoded to integer chromosome

      integer :: i
      integer :: j
      integer :: k
      integer :: aux1
      integer :: aux2

      z = 0
      do i = 1, ps
         aux1 = 1
         do j = 1, nu
            aux2 = lu(j) - 1
            do k = aux1, aux1+lu(j)-1
               z(i,j) = z(i,j) + pop(i,k)*2**aux2
               aux2 = aux2 - 1
            end do
            aux1 = aux1 + lu(j)
         end do
      end do

   end subroutine get_integer_value

   !============================================================================

   subroutine selection_roullete(ps, fit, sumfit, mate) ! last parameter is output
      implicit none
      integer, intent(out) :: mate     !< number of the individual for crossover
      integer, intent(in) :: ps        !< population size
      real(8), intent(in) :: sumfit    !< sum of fitness
      real(8), intent(in) :: fit(ps)   !< fitness of the individuals

      integer :: i
      real(8) :: rnd
      real(8) :: sumaux
      real(8) :: aux

      call random_number(rnd)
      aux = rnd*sumfit

      sumaux = 0.d0
      do i = 1, ps
         sumaux = sumaux + fit(i)
         if (sumaux > aux) exit
      end do

      mate = i

   end subroutine selection_roullete

   !============================================================================

   subroutine selection_roullete2(ps, pmin, fit, mate) ! last parameter is output
      implicit none
      integer, intent(out) :: mate     !< number of the individual for crossover
      integer, intent(in) :: ps        !< population size
      real(8), intent(in) :: pmin      !< minimum probability in crossover
      real(8), intent(in) :: fit(ps)   !< fitness of the individuals

      integer :: i
      real(8) :: rnd
      real(8) :: sumaux
      real(8) :: minfit    !< lower fitness
      real(8) :: maxfit    !< biggest fitness
      real(8) :: aux1
      real(8) :: faux(ps)

      minfit = minval(fit)
      maxfit = maxval(fit)

      call random_number(rnd)
      do i = 1, ps
         faux(i) = (1 - pmin)*(fit(i) - minfit)/(maxfit - minfit) + pmin
      end do

      aux1 = rnd*sum(faux)

      sumaux = 0
      do i = 1, ps
         sumaux = sumaux + faux(i)
         if (sumaux > aux1) exit
      end do

      mate = i

   end subroutine selection_roullete2

   !============================================================================

   !> \brief Pass chromosomes to the temporary population with performing crossover
   subroutine crossover_simple(mate1, mate2, ni, ps, tl, pop, newpop) ! last parameter is output
      implicit none
      integer, intent(in) :: mate1  !< number of the individual for crossover
      integer, intent(in) :: mate2  !< number of the individual for crossover
      integer, intent(in) :: ni     !< number of the individual
      integer, intent(in) :: ps     !< population size
      integer, intent(in) :: tl     !< total lenght of unknowns
      integer, intent(in) :: pop(ps,tl)     !< population of chromosomes
      integer, intent(out) :: newpop(ps,tl) !< new population of chromosomes

      integer :: j         !< counter
      integer :: crosssite !< number of genes crossed
      real(8) :: rnd       !< random number

      call random_number(rnd)

      crosssite = int((tl - 1)*rnd) + 1

      do j = 1, crosssite
         newpop(ni,j) = pop(mate1,j)
         newpop(ni+1,j) = pop(mate2,j)
      end do

      do j = crosssite+1, tl
         newpop(ni,j) = pop(mate2,j)
         newpop(ni+1,j) = pop(mate1,j)
      end do

   end subroutine crossover_simple

   !============================================================================

   !> \brief Pass chromosomes to the temporary population without performing crossover
   subroutine nocrossover(mate1, mate2, ni, ps, tl, pop, newpop) ! last parameter is output
      implicit none
      integer, intent(in) :: mate1  !< number of the individual for crossover
      integer, intent(in) :: mate2  !< number of the individual for crossover
      integer, intent(in) :: ni     !< number of the individual
      integer, intent(in) :: ps     !< population size
      integer, intent(in) :: tl     !< total lenght of unknowns
      integer, intent(in) :: pop(ps,tl)     !< population of chromosomes
      integer, intent(out) :: newpop(ps,tl) !< new population of chromosomes

      integer :: j

      do j = 1, tl
         newpop(ni,j) = pop(mate1,j)
         newpop(ni+1,j) = pop(mate2,j)
      end do

   end subroutine nocrossover

   !============================================================================

   !> \brief Mutation over the temporary population
   subroutine mutation(ps, tl, pm, newpop) ! last parameter is output
      implicit none
      integer, intent(in) :: ps              !< population size
      integer, intent(in) :: tl              !< total lenght of unknowns
      real(8), intent(in) :: pm              !< mutation probability
      integer, intent(out) :: newpop(ps,tl)  !< new population of chromosomes

      integer :: i
      integer :: j
      real(8) :: rnd

      do i = 1, ps
         do j = 1, tl

            call random_number(rnd)

            if (rnd <= pm) then
               if (newpop(i,j) == 1) then
                  newpop(i,j) = 0
               else
                  newpop(i,j) = 1
               end if
            end if

         end do
      end do

   end subroutine mutation

   !============================================================================

   !> \brief Replace the old population completely by the new one
   subroutine replace(ps, tl, newpop, pop) ! last parameter is output
      implicit none
      integer, intent(in) :: ps              !< population size
      integer, intent(in) :: tl              !< total lenght of unknowns
      integer, intent(in) :: newpop(ps,tl)   !< new population of chromosomes
      integer, intent(out) :: pop(ps,tl)     !< population of chromosomes

      integer :: i, j

      do i = 1, ps
         do j = 1, tl
            pop(i,j) = newpop(i,j)
         end do
      end do

   end subroutine replace

   !============================================================================

   subroutine get_fittest(ps, tl, fit, pop, newpop) ! last parameter is output
      implicit none
      integer, intent(in) :: ps              !< population size
      integer, intent(in) :: tl              !< total lenght of unknowns
      integer, intent(in) :: pop(ps,tl)      !< population of chromosomes
      integer, intent(out) :: newpop(ps,tl)  !< new population of chromosomes
      real(8), intent(in) :: fit(ps)         !< fitness of chromosome

      integer :: i
      integer :: ibest     !< index of the fittest
      real(8) :: fittest   !< biggest fitness

      ibest = 1
      fittest = fit(ibest)
      do i = 1, ps
         if (fit(i) > fittest) then
            fittest = fit(i)
            ibest = i
         end if
      end do

      newpop(1,:) = pop(ibest,:)

   end subroutine get_fittest

   !============================================================================

   subroutine convert_int_to_char3(nmbr, char3) ! lasta parameter is output
      implicit none
      character*3, intent(out) :: char3
      integer, intent(in) :: nmbr

      if (nmbr < 10) then
         write(char3,"('00',i1)") nmbr
      else if (nmbr < 100) then
         write(char3,"('0',i2)") nmbr
      else
         write(char3,"(i3)") nmbr
      end if

   end subroutine convert_int_to_char3

   !============================================================================

   subroutine plot_statistics(sname, folderout)
      implicit none
      character(100), intent(in) :: sname     !< simulations name
      character(100), intent(in) :: folderout !< folder the for output files

      open(48, file = trim(folderout) // trim(sname) // "-statistics.plt")

      write(48,*) "set term postscript color"
      write(48,*) "set out '", trim(folderout), trim(sname), "-statistics.eps'"
      write(48,*) "set style data linespoints"
      write(48,*) "set time"
      write(48,*) "set title '", trim(sname), "'"
      write(48,*) "set xlabel 'generation  (g)'"
      write(48,*) "set ylabel 'fitness'"
      write(48,*) "set grid"
      write(48,*) "plot '", trim(folderout), trim(sname), &
         "-statistics.txt' using 1:3 pt 4 lc 1 title 'Fittest',\"
      write(48,*) "     '", trim(folderout), trim(sname), &
         "-statistics.txt' using 1:2 pt 6 lc 3 title 'Mean'"
      write(48,*) "exit"

      close(48)

      call system("gnuplot -persist "// trim(folderout) // trim(sname) // &
         "-statistics.plt")

   end subroutine plot_statistics

   !============================================================================

   subroutine save_backup(folderout, sname, tcpuo, g, ffit, ps, ng, elitism, &
         nu, tl, pc, pm, pmin, lu, pop, rmin, rmax)
      implicit none
      character(100), intent(in) :: sname     !< simulations name
      character(100), intent(in) :: folderout !< folder the for output files
      character(100), intent(in) :: ffit      !< name of executable for fitness calculation

      integer, intent(in) :: g            !< generation
      integer, intent(in) :: ng           !< number of generations
      integer, intent(in) :: ps           !< population size
      integer, intent(in) :: nu           !< number of unknowns
      integer, intent(in) :: tl           !< total lenght of unknowns
      integer, intent(in) :: elitism      !< 0=off; 1=on
      integer, intent(in) :: lu(nu)       !< lenght of unknowns
      integer, intent(in) :: pop(ps,tl)   !< population of chromosomes

      real(8), intent(in) :: tcpuo        !< Accumulated CPU time
      real(8), intent(in) :: pc           !< crossover probability
      real(8), intent(in) :: pm           !< mutation probability
      real(8), intent(in) :: pmin         !< minimum probability in crossover
      real(8), intent(in) :: rmin(nu)     !< minimum real value of the unknowns
      real(8), intent(in) :: rmax(nu)     !< maximum real value of the unknowns

      integer :: i

      open(23, file = trim(folderout) // trim(sname) // "-backup.txt")

      if (len(trim(sname)) <= 15) then
         write(23,"(a15, a)") trim(sname), " = sname:   Simulation name"
      else
         write(23,"(a, a)") trim(sname), " = sname: Simulation name"
      end if

      write(23,"(1pe15.8, a)") tcpuo, " = tcpuo:   Accumulated CPU time"

      if (len(trim(ffit)) <= 15) then
         write(23,"(a15, a)") trim(ffit), &
            " = ffit:    Name of executable for fitness calculation"
      else
         write(23,"(a, a)") trim(ffit), &
            " = ffit: Name of executable for fitness calculation"
      end if

      write(23,"(i15, a)") ps, " = ps:      Population size (even number)"
      write(23,"(i15, a)") ng, " = ng:      Number of generations"
      write(23,"(1pe15.8, a)") pc, " = pc:      Crossover probability [0:1]"
      write(23,"(1pe15.8, a)") pm, &
         " = pm:      Mutation probability / (lu(1) + ... + lu(nu)) [0:1]"
      write(23,"(1pe15.8, a)") pmin, &
         " = pmin:    Minimum probability in crossover [0:1]"
      write(23,"(i15, a)") elitism, " = elitism: (0=off; 1=on)"
      write(23,"(i15, a)") nu, " = nu:      Number of unknowns"

      do i = 1, nu
         write(23,"(i15, a, i1, a, i1, a)") lu(i), &
            " = lu(", i, "):    Length of the ", i, " unknown"
         write(23,"(1pe15.8, a, i1, a, i1, a)") rmin(i), &
            " = rmin(", i, "): Lower value for the ", i, " unknown"
         write(23,"(1pe15.8, a, i1, a, i1, a)") rmax(i), &
            " = rmax(", i, "): Biggest value for the ", i, " unknown"
      end do

      write(23,"(i15, a)") g, " = g:       Last generation"

      write(23,*)
      do i = 1, ps
         write(23,"(1000(i2))") pop(i,:)
      end do
      write(23,*)

      close(23)

   end subroutine save_backup

   !============================================================================

   subroutine convert_real_to_time_format(rtime, ftime)
      implicit none
      character*10, intent(out) :: ftime  !< time formatted
      real*8, intent(in) :: rtime         !< time in seconds

      character*2 :: char2
      integer :: hora
      integer :: min
      integer :: seg

      hora = int(rtime/3600)
      min = int((rtime - hora*3600)/60)
      seg = int(rtime - hora*3600 - min*60)

      write(ftime, "(i4)") hora

      write(char2, "(i2)") min
      if (min < 10) char2 = "0" // adjustl(char2)
      ftime = trim(ftime) // ":" // char2

      write(char2, "(i2)") seg
      if (seg < 10) char2 = "0" // adjustl(char2)
      ftime = trim(ftime) // ":" // char2

   end subroutine convert_real_to_time_format

   !============================================================================

   subroutine write_results(tcpu, nu, ps, fit, r)
      implicit none
      integer, intent(in) :: nu        !< number of unknowns
      integer, intent(in) :: ps        !< population size
      real(8), intent(in) :: tcpu      !< total CPU time
      real(8), intent(in) :: fit(ps)   !< fitness of the population
      real(8), intent(in) :: r(ps,nu)  !< decoded to real unknowns

      character(10) tcpuf
      integer :: j

      call convert_real_to_time_format(tcpu, tcpuf)

      write(*,*)
      write(*,*) "SOLUCAO:"
      do j = 1, nu
         write(*,"('R(', i1, ') = ', 1pe15.8)") j, r(1,j)
      end do

      write(*,*)
      write(*,"('FITNESS: ', 1pe15.8)") fit(1)
      write(*,*)

      write(22,*)
      write(22,*) " ==================  SOLUTION  ================== "
      write(22,*)

      do j = 1, nu
         write(22,"(1pe15.8, a, i1, a, i1, a)") r(1,j), &
            " = r(", j, "):    The best value for the ", j, " unknown"
      end do

      write(22,*)
      write(22,"(1pe15.8, a)") fit(1), " = fittest: The biggest fitness found"
      write(22,*)

      write(22,"(1pe15.8, a)") tcpu, " = tcpu:    Total CPU time [s]"
      write(22,"(a15, a)") tcpuf, &
         " = tcpuf:   Total CPU time formatted [hh:mm:ss]"
      write(22,*)

      close(22)

   end subroutine write_results

   !============================================================================

   subroutine plot_history(sname, folderout, nu, ps, ng, xmin, xmax)
      implicit none
      character(100), intent(in) :: sname       !< simulations name
      character(100), intent(in) :: folderout   !< folder the for output files
      integer, intent(in) :: nu        !< number of unknowns
      integer, intent(in) :: ps        !< population size
      integer, intent(in) :: ng        !< maximal number of generations
      real(8), intent(in) :: xmin(nu)  !< low boundary constraints
      real(8), intent(in) :: xmax(nu)  !< high boundary constraints

      integer :: i
      integer :: j
      integer :: k
      integer :: m
      real(8) :: hx
      real(8) :: hy
      real(8) :: xx
      real(8) :: yy
      real(8) :: par1(ps*ng)
      real(8) :: par2(ps*ng)
      real(8) :: fitness(ps*ng)

      if (nu == 2) then

         open(46, file = trim(folderout) // trim(sname) // "-history.txt")

         read(46,*)
         do i = 1, ps*ng
            read(46,*) j, k, fitness(i), par1(i), par2(i)
         end do

         close(46)

         open(45, file = trim(folderout) // trim(sname) // "-field.txt")

         write(45,*) "#   parameter 1      parameter 2          fitness"

         m = 256
         hx = (xmax(1) - xmin(1))/(m - 1)
         hy = (xmax(2) - xmin(2))/(m - 1)

         do j = 1, m

            do i = 1, m

               xx = xmin(1) + (i - 1)*hx
               yy = xmin(2) + (j - 1)*hy

               call get_nearest(xx, yy, hx, hy, ps*ng, par1, par2, k)

               write(45,"(3(1pe15.8, 2x))") xx, yy, fitness(k)

            end do

            write(45,*)

         end do

         write(45,*)

         close(45)

      end if

      open(47, file = trim(folderout) // trim(sname) // "-history.plt")

      write(47,*) "set term postscript color"
      write(47,*) "set out '", trim(folderout), trim(sname), "-history.eps'"
      write(47,*) "set time"
      write(47,*) "set title 'history of ", trim(sname), "'"

      if (nu == 1) then

         write(47,*) "set style data points"
         write(47,*) "set grid"
         write(47,*) "set xrange [", xmin(1), ":", xmax(1), "]"
         write(47,*) "set xlabel 'parameter'"
         write(47,*) "set ylabel 'fitness'"
         write(47,*) "plot '", trim(folderout), trim(sname), &
            "-history.txt' using 4:3 pt 4 lc 1 title ''"

      else if (nu == 2) then

         write(47,*) "set style data pm3d"
         write(47,*) "set ticslevel 0"
         write(47,*) "set view map"
         write(47,*) "set pm3d map"
         write(47,*) "set palette"
         write(47,*) "set xrange [", xmin(1), ":", xmax(1), "]"
         write(47,*) "set yrange [", xmin(2), ":", xmax(2), "]"
         write(47,*) "set xlabel 'parameter 1'"
         write(47,*) "set ylabel 'parameter 2'"
         write(47,*) "set zlabel 'fitness'"
         write(47,*) "splot '", trim(folderout), trim(sname), &
            "-field.txt' using 1:2:3 title ''"

      end if

      write(47,*) "exit"

      close(47)

      if (nu == 1 .or. nu == 2) then
         call system("gnuplot -persist "// trim(folderout) // trim(sname) // &
            "-history.plt")
      end if

   end subroutine plot_history

   !============================================================================

   subroutine get_nearest(xx, yy, hx, hy, m, par1, par2, k)
      implicit none
      integer, intent(in) :: m
      integer, intent(out) :: k
      real(8), intent(in) :: xx
      real(8), intent(in) :: yy
      real(8), intent(in) :: hx
      real(8), intent(in) :: hy
      real(8), intent(in) :: par1(m)
      real(8), intent(in) :: par2(m)

      integer :: i
      real(8) :: ld
      real(8) :: dist

      k = 1
      ld = sqrt(((xx - par1(1))/hx)**2 + ((yy - par2(1))/hy)**2)

      do i = 1, m

         dist = sqrt(((xx - par1(i))/hx)**2 + ((yy - par2(i))/hy)**2)

         if (dist < ld) then
            k = i
            ld = dist
         end if

      end do

   end subroutine get_nearest

   !============================================================================

end program gapp
